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PACS 47 . 57 . Gc - Granular flow 

PACS 81 .05.Rm - Porous materials; granular materials 
PACS 47 . 50 . -d - Non-Newtonian fluid flows 

Abstract. - We discuss the advantages and results of using a mixing- length, compressible model 
to account for shear banding behaviour in granular flow. We formulate a general approach based 
on two function of the solid fraction to be determined. Studying the vertical chute flow, we show 
that shear band thickness is always independent from flowrate in the quasistatic limit, for Coulomb 
wall boundary conditions. The effect of bin width is addressed using the functions developed by 
Pouliquen and coworkers, predicting a linear dependence of shear band thickness by channel width, 
while literature reports contrasting data. We also discuss the influence of wall roughness on shear 
bands. Through a Coulomb wall friction criterion we show that our model correctly predicts the 
effect of increasing wall roughness on the thickness of shear bands. Then a simple mixing-length 
approach to steady granular flows can be useful and representative of a number of original features 
of granular flow. 
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Introduction. — Granular materials are ubiquitous 
in everyday's life, but a satisfying comprehension of their 
complex dynamics has not yet been achieved, despite the 
efforts of condensed matter science [1] . From the latin poet 
Lucretius to nowadays, we have experienced a perspective 
reversal; while, in his Dc Rerum Natura, he used seeds as 
, an example to explain the "flowability" of liquids, today 
'we try and use hydrodynamic analogies (e.g. [2], [3], [4], 
[5]) to explain the flow of granular materials. However, 
towards those poppy seeds our various theories remain al- 
most a divinatory exercise. 

Many of the theoretical approaches appeared in last de- 
cade literature neglect the compressibility of granular ma- 
terials, assuming it as an incompressible fluid with p 
const. From a phenomenological point of view, we believe 
that dilatancy is a requirement for shearing a granular ma- 
terial, in other words, the material has to dilate in order 
to shear. Accordingly, we expect that neglect of dilatancy 
loses an important part of the granular flow physics. 
In this work we formulate a simple model explicitly in- 
volving the solid fraction influence on the flow properties 
of granular materials. The model aim at being a gener- 
alization of the one by Pouliquen et al. [6] [7] based on 
the dimensionless parameter /, considered as the ratio be- 
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tween shearing time and rearranging time due to pressure. 
In Pouliquen's formulation, the solid fraction is derived 
from / as being linearly depending on it. Here we reverse 
the formulation assuming the solid fraction <j) to be the 
critical variable, instead of /, to reestablish the physical 
relevance of the dilation of the medium to determine the 
flow features. Giving the model an appropriate account of 
the solid fraction can become important for those geome- 
tries (like silos) in which <f> varies significantly (more than 
10%) all over the flow section. In silos flow can be seen 
to originate from fluidization due to the injection of voids 
from the exit hole, where solid fraction is quite different 
from its value in the core. 

In this perspective, we illustrate the 0-based model deriva- 
tion and its application to the vertical chute arrangement, 
to verify the constitutive relations proposed. 

The model. — The model outlined here is formulated 
for 2D, steady granular flows. The relevant equations are 
momentum balance with its two components, and the 
equation of continuity. 

As a fundamental assumption, we consider the flow 
structure to be solvable with the steady, compressible 
Navier-Stokes (N-S) equations, with a solid fraction, 
pressure and shear-dependent viscosity, n addition, we 
need a constitutive equation for ry, and we saturate the 
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degree of freedom introduced by <fi with a sort of Equation 
of State that involve pressure. 
Accordingly, equations are: 



V • pu = 



(1) 



V • (puu) 



V • (pvu) 
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du dv 
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dp d ( dv\ d f du dv\ 

(3) 

where p is the local density of the compressible pseudo- 
homogeneous medium, which indeed is a two phase mix- 
ture: 

P = Pp4> + Pint(l ~ 4>) (4) 

being p p and pi nt material and interstitial fluid densities, 
respectively. If the interstitial fluid is air, given that 
Pair « P P , we can approximate 



P P <i 
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Constitutive relations. We assume a rheological law 
based on dimensional analysis, like Prandtl's approach to 
turbulent flows, as discussed by Ertas and Halscy [8]. Be- 
cause of the eminent precursor, we will call it "mixing 
length approach" . 

In this perspective the apparent viscosity of the medium 
is locally formulated as: 
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where the unique timescale is pi ; L is a characteristic 
length, that has to be function of d and <f> only, with a 
generic relation of the form: 
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The function /(</>) is not known so far, but some features of 
it may be prescribed: it should diverge when 4> — > 4> ma x, 
to limit the material flow (that becomes 'jammed'), i.e. 
r\ —* 00. To achieve this limit r\ should diverge faster than 
|7| _1 , as it can be easily seen from Eq. 6. Interestingly, 
for f((f>) — 1 Eq. 7 reduce to Bagnold scaling for shear 
stress (valid for rapid granular flows), providing a further 
requirement that f(<p — > 0) = 1. However, the present 
work addresses dense flow of granular material, and we 
are not interested, at the moment, in the liquid-gas like 
transition. 

We also need a relation between pressure and solid frac- 
tion, which is similar to an Equation of State (EoS). As- 
suming shear rate plays the role of temperature in a gas, 
and <j) acts through a geometrical (excluded volume) func- 
tion h{4>) to be specified, we can use dimensional analysis 
the obtain the following EoS: 



p = Pp h{4>) (| 7 | df 



(8) 



To keep pressure finite when shear rate vanishes, h{cj)) has 
to diverge when <f> — ► 4> max (although value and physical 
meaning of (f> max is still a matter of debate [9]). 
Dimensional analysis is broadly used in granular flow mod- 
eling attempts, starting from Bagnold's works [10]. The 
formulation of Josserand et al. [2] [9] uses dimensional 
analysis with Coulomb friction to develop a constitutive 
relation for shear stress that is composed by a rate de- 
pendent part and a rate independent one, and where the 
isotropic part of the stress tensor is related to solid fraction 
by means of entropic considerations. We express normal 
and shear stresses according to Pouliquen [6] [7] , with the 
difference that we use explicitly solid fraction as the key 
variable, instead of dimensionless shear rate. Note that 
these laws, are very similar to those developed from hy- 
drodynamic analogies [4] [5], where granular temperature 
is used to represent the local mobility of the medium. We 
prefer our simple closure, based on 4> and an EoS for it, 
also because granular temperature is a variable which is 
difficult to measure and then correctly validate. 
Rearranging Eq. 6-7, we obtain: 
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For sake of simplicity, we replaced the ratio f/h with G{<p) 
and introduce: 



F(4>) = [h(cj>)] 



-1/2 



(10) 



as a simple replacement, provided h always appears in this 
form in the following developments of NS eqs. Note that 
G must vanish if rj-j = r — > 0. 

Applying the model to the vertical chute. 

We choose the vertical chute configuration as a stan- 
dard benchmark for model evaluation. Original flow struc- 
tures, principally related to the width of shear zones, can 
be found in the chute flow, like in Couette cells. Refer- 
ence for these configurations is the well known paper by 
GDR MiDi [11]. Broadly speaking, the material flows in 
a plug-like fashion in the central part of the chute, while 
it is sheared near the wall. The extent of shear bands ap- 
parently approaches a typical dimension, of order 10-15 
particle diameters. Predicting shear bands' thickness is a 
benchmark for all models applied to the chute and Couette 
flow [7]. 

Vertical chute equations. A scheme of the chute is 
given in Fig. 1. For the steady vertical chute 2D flow, 
N-S equations simplify thanks to: 
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Figure 1: Vertical chute scheme. 

We also included Janssen effect, assuming pressure does 
not vary in the vertical direction. Eq. 12 additionally 
states that pressure will not vary horizontally either. The 
other equation can be integrated. With y originating from 
the center of symmetry and directed to the wall, as shown 
in Fig. 1, we have: 



du 
dy 



(13) 



to be replaced in the s-momentum balance in combination 
with Eq. 9, to give: 

d ( d u\ 

d~V V [d^J + Pp(t>9 = ~ V d~y + Pp4>9 = ° (W) 



or 



G(4>) = — 
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G(y = 0) = follows from Eq. 9 and symmetry, which 
requires the shear stress at the centerline to vanish. 
From Eq. 15 we expect to identify <j>(y) provided p and 
an invertible form of G are given. At the same time, u(y) 
can be obtained from the EoS: 
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u(y) = u sHp + f F{4>) dy (17) 



or 



where b is the half-width of the channel. So far, the un- 
known functions <j>(y) and u(y) can be formally obtained 

as: 

J G(0) = effiWv 

IvTpv rb 



u(y) 



u slip + ^f°F(ct>)dy 



(18) 



but in practice F(<p) and G{4>) must be specified, and also 
the pressure calculated. 

The continuity equation can be used in its integral form, 
correlating local profiles to the total mass flowrate, M. 
Experiments are easily performed with constant flowrate, 
either controlled by a simple hole in the bottom of the silo, 
or using a moving plate with fixed velocity. Accordingly: 
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Developing, we can formulate an expression for the slip 
velocity by using Eq. 17: 
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On boundary conditions. One of the most critical issue 
in granular flow simulation is identification and applica- 
tion of boundary conditions. For the solid fraction, we 
will assume that in the central zone of the chute material 
reaches random close packing, or some other critical value 
of (p that leads pressure to diverge [9]. 
One constraint on velocity can be formulated as an integral 
condition, by fixing flowrate as done above. In addition, 
we must speculate on the interaction between the granular 
assembly as a continuum and the walls. The simplest view 
used in literature is assuming a layer of particles glued at 
the walls, for which we can use a no-slip boundary condi- 
tion. This assumption is attracting for its simplicity but 
we believe it requires caution in its application. We are 
skeptic that the continuum averaged interaction between 
nearest particles is the same in the bulk and in the layer of 
particles facing the glued ones. In this perspective, exper- 
imental investigation and critical theoretical speculations 
have to be done. A viable alternative to no-slip assump- 
tion is the Coulomb criterion at the wall: 



,,tanS 



(22) 



where 8 is a characteristic wall friction angle. In case of 
particle artificially fixed at the wall, this means assuming 
them as a wall, with a specific roughness measurable by 
its own 8. Combining Coulombs law with Eq. 9 we obtain 
a condition on the solid fraction: 



|7|, 



-G(<j> w )\j\ w =pG(<j> w ) (23) 



which leads to 

G{(j> w ) = tanS (24) 

given that o~ w = p. Very important, with Coulombs crite- 
rion the slip velocity is not zero, and has to be determined 
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from Eq. 18 using flowrate. In addition Eq. 24 allows to 
calculate the pressure (normal stress), provided that: 



Gw = G (y = b)=W f 
P Jo 



which gives: 



P = 



ave y v 



tanS 



tanS 



(25) 



(26) 



where the average solid fraction, defined by J Q 4>dy = 
b(j) aV e, has been introduced. 

Combining Eqs. 26, 21, and 18, the velocity profile can 
be explicitly written as: 



u(y) = 



M 
2p p 



>/<f>avegb 



dvtanS Jo 
dytanS 



<j)dy 



F(<P)dyj dy 



F(^)dy 



The result of Eq. 26 together with Eq. 15 leads to: 
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stating that in the approximation of 
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(28) 



(29) 



or G is a linear function of y, which provides a consistence 
criterion for the identification of the unknown function G. 

Interestingly, our model, based also on Coulomb wall 
criterion, predicts the invariance of the velocity profiles 
with flowrate in the quasistatic limit (where <fi w <j) a ve is 
valid). In other words, the scaled velocity profile: 



u = 



J b y F{ct>(y))dy 



(30) 



u m ax u s iip F((p(y))dy 

does not depend on flowrate, that influences only the slip 
tangential velocity. In this limit, the solid fraction profile 
also does not depend on flowrate. We underline that the 
result is independent of the particular formulation of the 
functions F and G. Also with Pouliquen's formulation for 
the effective friction coefficient (which can be seen as a 
particular choice for F and G), but with a Coulomb slip 
criterion, the independency of shear bands from flowrate 
is obtained. This is indeed a result supporting the mixing 
length approach. Profiles arc flowrate-indcpendent due 
to the nature of granular matter, which develops inter- 
nal stresses to sustain itself; in this perspective, taking 
into account stresses in the formulation of boundary con- 
ditions is necessary. However, far from the quasistatic 
limit, or conditions where 4> w <j) a ve is acceptable, flowrate 
can significantly affect <j>{y) and u(y) because of the close 
coupling of the two equations. 



Deriving expressions for F((f>) and G{4>). Following 
the work of Pouliquen and coworkers, the functions would 
take the form: 




Us + 



i /<t>+i 



where <f> is the scaled solid fraction given by: 



<t> 



4>r, 



(31) 



(32) 



The authors [7] acknowledged a major difficulty in the ap- 
plication to the vertical chute with no-slip at the walls; 
shear bands are not finite and of constant width in the 
quasistatic limit. We have already demonstrated that a 
Coulomb wall slip criterion can correct this. In the fol- 
lowing we will illustrate the results of our mixing length 
model including Coulomb wall slip criterion and F and G 
functions as in Eq. 31, for different chute width and wall 
roughness. 

Before that, note the analytical solution achievable in the 
quasistatic limit, obtained combining Eq. 29 and 31: 







>2" 



for y > n' s 
for y < (j,' s 



(33) 



where fi' are the corresponding parameters in Eq. 31 mul- 
tiplied by b/tanS. The scaled velocity profile, obtained 
combining Eq. 30,31 and 33, is indeed a simple function 
of y: 



A(y -b) + B ln(C - Dy) for y > M ' s 
1 for y< fi' s 



(34) 



where A,B,C,D are known constants, depending on model 
parameters, wall friction angle and channel width. 
In the following we perform an analysis based on these 
analytical results; far from the quasistatic limit, one can 
repeat the calculations using Eq. 28 instead of Eq. 29, in 
a numerical fashion. 

Shear bands thickness and chute width. We already 
demonstrated that the model predicts shear bands inde- 
pendent by flowrate. Here we explore the effect of bin 
width. 

Using the quasistatic assumption, wc calculate u from Eq. 
34, for different widths, b. Other model parameters must 
be given and they depend on the specific material choosen. 
for the purpose of illustrating the model, we used values 
determined by Jop et al. al. [6], collected in table 1. 
Fig. 2 shows that the model predicts a linear correlation 
between shear band width and channel width. The slope 
of the linear dependency may change with different ma- 
terials, but remains linear. It is frequently stated that 
the thickness of shear bands is expected to be indepen- 
dent from channel extension. However, literature reports 
data supporting (e.g. [12]) and contrasting [13] this state- 
ment. Our results agree with the experimental results by 
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Table 1: Parameters of the model [6]. 



h 


0.279 


(adim) 


Ms 


tan(20.9) 


(adim) 


M2 


tan(32.76) 


(adim) 



Nedderman and Lahoakul, but the issue requires further 
investigation of the vertical chute, in order to discriminate 
the applicability of a mixing-length model to this config- 
uration. 
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Figure 2: Non-dimensional axial velocity, solid fraction and 
shear bands thickness vs non-dimensional channel width, in 
chutes of different extension, (wall friction angle 8 = 25deg). 



Shear bands thickness and wall roughness. Wall 
roughness can influence the extension of the shear bands, 
according to our model. Adopting Coulomb criterion 
at the wall yields a simple expression of wall roughness, 
related to the wall friction angle S, while using a no-slip 
condition makes impossible to account for roughness 
within a continuum approach. Also in view of real scale 
application of a mixing-length model, wall friction must 
be correctly accounted for, and the wall friction angle 
is a wide-spread approach. Furthermore, real world 
applications aim at perfect wall slip, but are often in an 
intermediate situation, where wall roughness plays a role. 

Results according to Eq. 34 are given in Fig. 3, showing 
that the model predicts a dependence of shearing regions 
from wall roughness. However, the effect of S is approach- 
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Figure 3: Non-dimensional axial velocity, solid fraction and 
shear bands thickness vs non-dimensional channel width, in 
chutes with different wall friction (channel width = 50 particle 
diameters) . 



ing an asymptote. At small values slip occurs, whereas 
larger friction reduces its influence on the shearing bands. 
The enlargement shear zones with increasing wall rough- 
ness is supported also by DEM results of Prochnow [14]. 
Figure 4 illustrates qualitatively the comparison between 
DEM and our mixing-length, continuous model. The mix- 
ing length approach can capture the effect of increasing 
wall roughness predicted by DEM calculations by means 
of different sizes of particles fixed at the wall. 

Conclusions. — In this work we formulated and ap- 
plied a mixing-length, compressible continuum model, 
aiming at discussing shear band particular behaviour in 
granular flow. A general approach based on two functions 
of the solid fraction to be specified has been given. Ap- 
plication to the vertical chute problem demonstrates that 
for suitable wall boundary conditions (Coulomb friction), 
shear band thickness is independent of flowrate, in the 
quasistatic limit, whatever the form of the generic func- 
tions. The relevance of these boundary conditions in real 
applications was discussed, and a critical evaluation of ex- 
perimental practice was given. 

Further, we addressed the role of bin width, using expres- 
sion for the G((j)) and F(<p) functions based on reformula- 
tion of Literature results [7] to highlight the variable solid 
fraction. Results do not agree with the frequently assumed 
independence of shear band thickness from channel width, 
but agree with with Literature experimental data [13]. We 
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Figure 4: Rescaled velocity profiles from DEM simulations by 
Prochnow [14] and continuous, mixing-length model, (circles: 
D — 4 d, diamonds: D = d, squares: D — 0.5 d), while 
continuous model calculations use different <5 values. 

expect that additional experimental investigations might 
elucidate this controversy. Taking advantage of the inclu- 
sion of wall friction criterion instead of a no-slip boundary 
condition, we correctly predicted an asymptotic increase 
of shear bands extension with larger wall roughness. 
The analysis above proves that even a simplified mixing- 
length approach to steady granular flows can be useful 
and representative of a certain number of features, once 
the proper boundary conditions are used. 

* * * 

We acknowledge collaboration by Frangois Chevoir who 
kindly sent us his DEM results. 
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